环境设置

# 安装和加载 aplot
if (!require("aplot")) {
  if (!require("devtools")) install.packages("devtools")
  devtools::install_github("YuLab-SMU/aplot")
}

library(magrittr)
library(ggplot2)
library(ggrepel)
library(janitor)
library(aplot)
library(dplyr)
library(scales)

Sys.setenv(LANGUAGE = "en")
options(stringsAsFactors = FALSE)

数据加载和预处理

data <- read.csv("easy_input.csv", check.names = FALSE)
rownames(data) <- data[, 1]
colnames(data)[1] <- "Geneid"

data$neg_log10_p <- -log10(data$P.Value)
data$label <- ifelse(data$neg_log10_p > 60, rownames(data), "")

# 预计算全局范围,确保所有图使用相同的坐标系
x_range <- range(data$deltaBeta)
y_range <- c(0, 100)

head(data)
##                Geneid  P.Value  deltaBeta feature neg_log10_p label
## cg09403666 cg09403666 3.80e-11  0.5925933     IGR    10.42022      
## cg03433048 cg03433048 3.77e-11 -0.3084334  TSS200    10.42366      
## cg20636352 cg20636352 3.75e-11  0.5573535     IGR    10.42597      
## cg01854491 cg01854491 3.73e-11 -0.3970072   5'UTR    10.42829      
## cg25835179 cg25835179 3.70e-11 -0.3024118  TSS200    10.43180      
## cg26894523 cg26894523 3.69e-11  0.4178792    Body    10.43297

创建主火山图

create_volcano_plot <- function(df, x_lim, y_lim) {
  ggplot(df, aes(deltaBeta, neg_log10_p)) +
    geom_text_repel(aes(label = label), 
                    size = 3, max.overlaps = Inf,
                    box.padding = 0.3, point.padding = 0.3) +
    geom_point(aes(color = feature), alpha = 0.8, size = 1.2) +
    labs(x = "Methylation difference (beta-value)",
         y = bquote(~-log[10]~(italic("P-value")))) +
    theme_classic() +
    theme(
      axis.title = element_text(color = "black", size = 14, face = "bold"),
      axis.text = element_text(color = "black", size = 12),
      legend.position = "none",
      panel.border = element_rect(color = "black", fill = NA, size = 0.5),
      plot.margin = margin(5, 5, 5, 5)
    ) +
    scale_y_continuous(
      breaks = c(0, 25, 50, 75, 100),
      labels = c(0, 25, 50, 75, ""),
      limits = y_lim
    ) +
    scale_x_continuous(limits = x_lim) +
    coord_cartesian(xlim = x_lim, ylim = y_lim)
}

创建顶部百分比图

create_top_percentage_plot <- function(df, x_lim) {
  # 计算百分比
  percent_hyper <- sum(df$deltaBeta > 0) / nrow(df) * 100
  percent_hypo <- sum(df$deltaBeta < 0) / nrow(df) * 100
  
  # 创建百分比条数据 - 使用简单的数值坐标
  percentage_data <- data.frame(
    x = c(x_lim[1] + (0 - x_lim[1])/2, (x_lim[2] - 0)/2),
    y = c(0.5, 0.5),
    width = c(0 - x_lim[1], x_lim[2] - 0),
    height = c(1, 1),
    category = c("Hypomethylated", "Hypermethylated"),
    percentage = c(round(percent_hypo, 1), round(percent_hyper, 1)),
    fill_color = c("#BEBEBE", "#4D4D4D"),
    text_color = c("black", "white")
  )
  
  ggplot(percentage_data, aes(x = x, y = y)) +
    geom_tile(aes(width = width, height = height, fill = I(fill_color)), 
              color = "black", size = 0.5) +
    geom_text(aes(label = paste0(percentage, "%"), color = I(text_color)),
              fontface = "bold", size = 4) +
    theme_void() +
    theme(
      plot.margin = margin(0, 5, 0, 5),
      panel.background = element_blank()
    ) +
    scale_x_continuous(limits = x_lim) +
    scale_y_continuous(limits = c(0, 1)) +
    coord_cartesian(xlim = x_lim, ylim = c(0, 1))
}

创建右侧特征柱状图

create_feature_barplot <- function(df) {
  # 计算特征统计
  feature_stats <- df %>%
    count(feature) %>%
    mutate(percentage = round(n / sum(n) * 100, 1)) %>%
    arrange(desc(n)) %>%
    mutate(
      feature_clean = case_when(
        feature == "1stExon" ~ "1stExon",
        feature == "3'UTR" ~ "3'UTR", 
        feature == "5'UTR" ~ "5'UTR",
        feature == "Body" ~ "Body",
        feature == "IGR" ~ "IGR",
        feature == "TSS1500" ~ "TSS1500",
        feature == "TSS200" ~ "TSS200",
        TRUE ~ feature
      ),
      y_pos = row_number()
    )
  
  # 获取颜色映射
  colors <- scales::hue_pal()(length(unique(df$feature)))
  names(colors) <- sort(unique(df$feature))
  
  ggplot(feature_stats, aes(x = n, y = y_pos)) +
    geom_col(fill = "#BEBEBE", color = "black", width = 0.8) +
    geom_point(aes(color = feature), size = 3, x = 0) +
    geom_text(aes(label = paste0(percentage, "%")), 
              hjust = -0.1, fontface = "bold", size = 3) +
    geom_text(aes(x = -max(n)*0.1, label = feature_clean),
              hjust = 1, fontface = "bold", size = 3) +
    scale_color_manual(values = colors) +
    theme_void() +
    theme(
      legend.position = "none",
      plot.margin = margin(5, 5, 5, 0)
    ) +
    scale_x_continuous(limits = c(-max(feature_stats$n)*0.2, max(feature_stats$n)*1.3)) +
    scale_y_continuous(limits = c(0.5, nrow(feature_stats) + 0.5)) +
    coord_cartesian(
      xlim = c(-max(feature_stats$n)*0.2, max(feature_stats$n)*1.3),
      ylim = c(0.5, nrow(feature_stats) + 0.5)
    )
}

使用 aplot 组合图形

# 创建三个子图,使用统一的坐标限制
p_main <- create_volcano_plot(data, x_range, y_range)
p_top <- create_top_percentage_plot(data, x_range)  
p_right <- create_feature_barplot(data)

# 先测试单独的图是否正常
cat("测试主图:\n")
## 测试主图:
print(p_main)

cat("测试顶部图:\n")
## 测试顶部图:
print(p_top)

cat("测试右侧图:\n")
## 测试右侧图:
print(p_right)
## Warning: Removed 7 rows containing missing values or values outside the scale
## range (`geom_col()`).

# 使用 aplot 进行组合
cat("开始组合图形...\n")
## 开始组合图形...
# 方法1:分步组合
step1 <- insert_top(p_main, p_top, height = 0.15)
cat("第一步完成:添加顶部图\n")
## 第一步完成:添加顶部图
print(step1)

# 方法2:完整组合
volcano_combined <- insert_top(p_main, p_top, height = 0.15) %>%
  insert_right(p_right, width = 0.25)

cat("组合完成!\n")
## 组合完成!
print(volcano_combined)

保存图片

# 保存组合图
ggsave("multiVolcano_aplot_fixed.pdf", volcano_combined, 
       width = 12, height = 8, dpi = 300)

cat("图片保存完成:multiVolcano_aplot_fixed.pdf\n")
## 图片保存完成:multiVolcano_aplot_fixed.pdf

替代方案:基础包组合

# 如果上述方法仍有问题,使用基础包的组合方案
cat("尝试替代方案...\n")
## 尝试替代方案...
# 使用 grid 包进行图形组合
library(grid)
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
## 
##     combine
# 创建一个简单的组合图形
combined_plot <- grid.arrange(p_main, ncol = 1)

print(combined_plot)
## TableGrob (1 x 1) "arrange": 1 grobs
##   z     cells    name           grob
## 1 1 (1-1,1-1) arrange gtable[layout]

会话信息

sessionInfo()
## R version 4.5.1 (2025-06-13)
## Platform: x86_64-apple-darwin20
## Running under: macOS Monterey 12.7.6
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.5-x86_64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.5-x86_64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.1
## 
## locale:
## [1] zh_CN.UTF-8/zh_CN.UTF-8/zh_CN.UTF-8/C/zh_CN.UTF-8/zh_CN.UTF-8
## 
## time zone: Asia/Shanghai
## tzcode source: internal
## 
## attached base packages:
## [1] grid      stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] gridExtra_2.3   patchwork_1.3.2 tidyr_1.3.1     scales_1.4.0   
##  [5] dplyr_1.1.4     aplot_0.2.9     cowplot_1.2.0   janitor_2.2.1  
##  [9] ggrepel_0.9.6   ggplot2_4.0.0   magrittr_2.0.4 
## 
## loaded via a namespace (and not attached):
##  [1] yulab.utils_0.2.1  rappdirs_0.3.3     sass_0.4.10        generics_0.1.4    
##  [5] ggplotify_0.1.3    stringi_1.8.7      digest_0.6.37      evaluate_1.0.5    
##  [9] timechange_0.3.0   RColorBrewer_1.1-3 fastmap_1.2.0      jsonlite_2.0.0    
## [13] purrr_1.1.0        textshaping_1.0.4  jquerylib_0.1.4    cli_3.6.5         
## [17] rlang_1.1.6        crayon_1.5.3       withr_3.0.2        cachem_1.1.0      
## [21] yaml_2.3.10        tools_4.5.1        gridGraphics_0.5-1 vctrs_0.6.5       
## [25] R6_2.6.1           lifecycle_1.0.4    lubridate_1.9.4    snakecase_0.11.1  
## [29] stringr_1.5.2      ggfun_0.2.0        fs_1.6.6           ragg_1.5.0        
## [33] pkgconfig_2.0.3    pillar_1.11.1      bslib_0.9.0        gtable_0.3.6      
## [37] glue_1.8.0         Rcpp_1.1.0         systemfonts_1.3.1  xfun_0.53         
## [41] tibble_3.3.0       tidyselect_1.2.1   knitr_1.50         farver_2.1.2      
## [45] htmltools_0.5.8.1  rmarkdown_2.30     labeling_0.4.3     compiler_4.5.1    
## [49] S7_0.2.0